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Abstract. We consider a qubit subject to various independent control mechanisms and present a general strategy to identify 
both the internal Hamiltonian and the interaction Hamiltonian for each control mechanism, relying only on a single, fixed 
readout process such as d z measurements. 



1. INTRODUCTION 

Realizing the ultimate goal of quantum information processing, namely building a working quantum computer, is to a 
large extent a problem of finding ways to control the dynamics of a quantum system precisely. A crucial prerequisite 
for this task is one's ability to accurately determine of the dynamics of the physical system and its response to external 
(control) fields. Quantum process tomography (QPT), by providing a general procedure to identify the unitary (or 
completely positive) processes acting on a system, addresses this problem but does not solve it completely. 

One problematic aspect of QPT is the assumption that one can experimentally determine the expectation values 
of a complete set of observables, or at least perform arbitrary single qubit measurements on a register of n qubits. 
Most QIP proposals rely on a single readout process, i.e., measurement in a fixed basis. For example, qubits encoded 
in internal electronic states of trapped ions or neutral atoms are read out by quantum jump detection via a cycling 
transition; readout for solid-state qubits based on Cooper-pair boxes, Josephson junctions, or electrons in double-well 
potentials usually involves charge localization using single electron transistors or similar devices. Finally, solid state 
architectures based on electron or nuclear spin qubits are expected to be limited to a z measurements via spin-charge 
transfer. 

It is usually assumed that local projective measurements in a fixed basis are sufficient since arbitrary single qubit 
measurements can then be realized by performing a local unitary transformation before measuring to achieve a change 
of basis. However, implementing such a basis change requires precise knowledge of the dynamics of each individual 
qubit and its response to control fields in the first place, the very information we seek to determine experimentally, and 
which may be difficult to predict precisely based on theoretical models and computer simulations alone, in particular 
for systems that are sensitive to fabrication variance. A possible solution to this seemingly intractable problem was 
described in [1]. In the following we outline the general strategy, discuss various ways of extracting the system 
parameters from noisy experimental data, and illustrate the key steps using examples with simulated measurement 
data. 



2. GENERAL STRATEGY FOR HAMILTONIAN IDENTIFICATION 

The state of a two-level system can be mapped to a Bloch vector s, i.e., a real vector in R 3 , with pure states 
corresponding to points on the Bloch sphere, i.e., the surface of the unit ball. On timescales sufficiently short compared 
to the decoherence time, the evolution of the system is governed by a Hamiltonian, which can be written in terms of 
the Pauli matrices a* for * e {x,y,z}, 2H — d^I + d x a x + d y a y + d z a z , where do, d x , d y and d z are real constants. If the 
Hamiltonian remains constant for t q < t < fi, the Bloch vector undergoes a rotation about the axis d = (d x ,d y ,d z ). 




The length ||d|| of this vector determines rotation frequency; the rotation axis can be specified by a unit vector 
d = (sin0cos0,sin0sin0,cos0) r , i.e., by two angles 9 and (j) as shown in Fig. 1. To identify the parameters d x , 
d y and d z of the Hamiltonian it therefore suffices to determine the rotation frequency and the angles 9 and (p. Since do 
results only in an unobservable global phase factor, it can be ignored. 

If the system can be repeatedly initialized in a known state, e.g., one of the measurement basis states |0) or |1), 
and then measured after having evolved for progressively longer time periods, we can map the trajectory of the z- 
component of the Bloch vector, and extract the frequency and angle 9 of the rotation as shown in Fig. 1 (right). For 
some systems such as NMR-based schemes this may be sufficient as the phase relationship between d x and d y is fixed 
by the phase relationship between the control fields. In general, however, a second series of measurements is necessary 
to determine the horizontal angle (j> of the rotation axis with respect to a reference axis d, . This procedure is given in 
Ref. [1] and involves using the values of 9 and ||d|| determined in the first step to select a different initial state, and 
mapping its precession about the desired rotation axis. 

For a system subject to various control fields /,„ (e.g., fields associated with different control electrodes) in addition 
to its free evolution, we must determine both its internal Hamiltonian Ho and the interaction Hamiltonian H m for each 
independent control mechanism. Although we usually cannot determine the interaction Hamiltonians directly since 

(k) (k) 

we cannot switch off the internal dynamics, we can identify the rotation axis flo-i-m = ^o +/m "-m corresponding to 
the evolution of the system under the Hamiltonian Hq + fm^H m for a fixed control setting /„, = /jf . Repeating this 

(k) 

procedure for each available control field /„, with several control field settings /,„ then allows us to extract both the 
internal and interaction parts of the Hamiltonian provided that the dependence of the Hamiltonian on the control fields 
is linear. (Nonlinear field effects require additional correction terms, and we will exclude this case in this paper.) 



3. MAPPING z(t) AND EXTRACTING THE RELEVANT DATA 

A crucial factor in the Hamiltonian identification strategy outlined above is the mapping of the evolution of z(t) as 
the system precesses around a fixed rotation axis. The accuracy with which we can identify the relevant parameters 
such as the rotation frequency and angle 9 depends on the total length of time tf for which z{t) is mapped, the time 
resolution Af and the uncertainty of each data point z{tk) = (<7 z (f,t)}, which depends on the number of times N e each 
experiment is repeated to obtain the ensemble average (in addition to the frequency of measurement errors etc). The 
total number of measurements required to map z(t ) is thus Nj = N e tf/At. The choice of Af, tf and N e will depend both 
on the system to be characterized and the method of data extraction to be used. 

One possible approach is to use a small time step Af and a large number of repetitions N e to obtain a dense set of 
accurate data points for a period of time covering at least one quarter of the rotation period, and fit a cosine segment 
to the data. Since the time we have to monitor the evolution of z(t) is limited by the decoherence time of the system, 
this approach may be useful for systems that decohere rapidly because it requires only sampling over a short period of 
time. 
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Figure 2: The top graph to the left shows an 
example of noisy (simulated) measurement 

(2) 

data for rotations about do + f 2 &i with 
0.1. The dash-dot line indicates the 
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optimal sampling time tp determined by 
finding the maximum of the function P(tf) 
(bottom-left). The Fourier transform of the 
(truncated) data for < t < tf is shown in 
the bottom-right graph. Its zero-frequency 
component determines 9 — arccos ^/|F(0)|; 
and the location of the second peak gives the 

rotation frequency CO = ||do + &2\\- The 
estimates obtained, 9 = 0.80 and co = 0.28, 
are close to the actual values 9 = 0.7854 and 
(0 = 0.2828. 



Another approach, explained in detail in [1], is to use relatively coarse time sampling, with a moderate number 
of repetitions for each data point, over several rotation periods to obtain a rough estimate of the rotation frequency 
using the discrete Fourier transform, followed by a second step of acquiring more accurate additional data points in 
the region where the first minimum of z(t ) is expected based on the first estimate, and fitting a parabola to the new 
data to find the rotation frequency and declination of the rotation axis. 

A third alternative [2] is to eliminate curve-fitting altogether and extract all the required information directly from 
the Fourier transform. This method is rather elegant and does not require the high time resolution At and measurement 
repetition rates N e usually necessary for curve-fitting methods. All we need to avoid aliasing is that At be less than half 
the rotation period T. Ensuring that this condition is satisfied requires a rough a priori estimate of T but this should 
normally not be a problem. It also permits easy estimation of the accuracy of the parameters. However, since the 
frequency resolution Acq = 27r/f/, we must be able to map the data for a least two complete cycles to be able to extract 
the rotation frequency, and more cycles will be required to obtain a clearly defined peak in the frequency spectrum. 
Hence, this method will be most suited to systems whose decoherence time is sufficiently long to allow mapping of 
the evolution of z{t ) over several cycles. 



4. ILLUSTRATIVE EXAMPLE 

To demonstrate the procedure, we choose a test system with do = (0.2,0,0.1) r , di = (1, l,0) r and d2 = (0,0, l) r . 
Fig. 2 illustrates how we identify the rotation frequencies and declination angles for f\ = and fi =0.1 following the 
3rd approach outlined in Sec. 3. We sample z(t ) over several rotation periods with an intermediate time step (At = 0.25) 
and a small number of repetitions (N e = 10) for each measurement. We then obtain an estimate of the rotation frequency 
by taking the Fourier transform of the data and finding the frequency co p such that \F(co p )\ = max„>o |F(a)„)|, where 
F(co„) is the nth Fourier coefficient. Since the sampling period is usually not an integer multiple of the rotation period, 
the peak in the Fourier spectrum will tend to be unsharp, and our estimate inaccurate. To improve it, we compute 
the function P{tf) = [F tf (co p ) — F tf {co p -\) — F tf (co p+ \)\/[F tf {co p -i) + F tf (co p+ i)] where F tf is the Fourier transform 
of the truncated data for < t < tf and co p is the frequency where the first peak in the spectrum (excl. F(Q)) occurs, 
\F, f (co p )\ = max„>o \F tf (co„)\. P(tf) assumes a maximum when tf is an integer multiple of the rotation period, thus 
allowing us to find the optimal sampling time tf. Fig. 3 illustrates how we can find the horizontal angles 0, having 
identified the rotation frequencies and angles 9 for all control settings, using a local curve fitting approach similar to 
the second strategy outlined in Sec. 3 and described in Ref. [1]. Finally, Fig. 4 shows how we can extract do and d„, 
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by plotting the x, y and z-components of the rotation axes do + /,,, d,„ versus /,;, for m — 1 , 2 and fitting straight lines 
to the data. 



Figure 3: The system is initialized in state S\ = 
(cos/3,sin/3,0) r (here /3 = -1.0446) by rotating s = 
(0, 0, 1 ) T about the reference axis d,- (here do) by a suitable 

angle y/, and the precession of Si about do + f^di, whose 
frequency and 6 are already known, is mapped. Using a 
small number of data points (stars) we find the x-intercept 
of z(t) (circle), which allows us to estimate the location of 
the extrema of z(t ). Additional data points (dots) are then 
acquired in these regions, and parabolas fitted to the data. 
The desired angle <j> = — j5 — arcsin(ycos 8 / sin 9 r ), where 
7 = (Zmax ~ Zmin)/2 and 8 = K - (a min + a max )/2 are de- 
termined by the vertices of the parabolas, and the values /3 
and 9 r (known from part 1). We obtain = 0.34, which is 
close to the actual value <j> — 0.3218. 
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Figure 4: Having determined the rotation frequen- 
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cies CBfa. m and the angles 9^ m and <p^ m of the 
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rotation axes do + fm d,„ for various values of the 
controls f\ and fi, we convert the data into Carte- 
sian coordinates, plot the values of the x, y and z- 
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components of the axes do + f,„ d m for m— 1,2, re- 
spectively, and fit straight lines. The y-intercepts de- 
termine the x, y and z-components of do, the slopes 
those of d m , for m — 1,2. For the data shown, we 
have 

do" = (0.1986, 0.0048, 0.0979) r 
df = (0.9884, 1.0163, 0.0087) r 
df = (0.053 1,0.0246, 0.98 19) r 



(*) 



The distances ||df,f -d^'|| of 0.0054, 0.0218 and 
u 5 0.0613, respectively, compare favorably to the 3% 
readout error rate of the simulated experiments. 
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